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I. INTRODUCTION 

We will consider the construction of the absorbing 
boundary conditions for the Schrodinger type equations 
(subscripts are used for the derivative with respect to the 
corresponding variable) 

iux + (3{x)uyy + i^{x,y)u ^ , (1) 

which arise in numerous quantum-mechanical problems 
and as approximate models for the wave propagation in 
the parabolic equation method and its generalizations [ij . 
For calculating solutions of this equation in unbounded 
domains it is essential to introduce boundary conditions 
at the boundaries of the computational domain which 
model the free transmission of the waves through these 
boundaries. As these boundary conditions must mini- 
mize the amplitudes of waves reflected from boundaries, 
they are called absorbing boundary conditions [2| . 

There is a significant number of works where the prob- 
lem of constructing of such boundary conditions was 
considered. The main approaches in these works con- 
sist either in factorization of the differential operator of 
the equation under consideration into pseudodifferential 
factors, each of which describes the unidirectional wave 
propagation 0, 13j Q and use the differential approxima- 
tions of these factors for the formulation of the boundary 
conditions, or in formulation of the absorbing boundary 
conditions as the matching condition with the free space 
solution outside of the computational domain (for the 
Eq. ([T]) see the paper 

As is known, the approximate description of the uni- 
directional propagation of the waves can be obtained by 
the generalized multiple-scale method 0), which in par- 
ticular cases gives the same results as the WKB or ray 
method. In this approach the algebraic factorization of 
the Hamilton-Jacobi equation replaces the factorization 
of the differential operator. In the simplest case this ap- 
proach was reported in 



II. DERIVATION OF THE BOUNDARY 
CONDITIONS 

To apply the multiple-scale method [7] to the Eq. ([T|) 
we introduce the slow variables X = ex, Y = ey, the fast 
variable rj = {l/e)6{X, Y), where e is a small parameter, 
change the partial derivatives in Eq. ([T|) for the prolonged 
ones by the rules 

d d ^ d d d ^ d 

dx dX ^ dr]' dy dY ^ dr] ' 

and substitute in the obtained equation the expansion 

u = uq + eul -|- . . . . 

Equating coefficients of like powers of e, we obtain at 
0(e°) the representation 

uo = A^{X,Yy'^ + B,,{X,Y)e-'\ 

and the Hamilton-Jacobi equation for the phase function 

e 

9x + P{BYf -v = Q. (2) 
Later we will consider the one-way part of this solution 

Mo = Aoe"' = Ao{X, Y) exp(i e{X, Y)) . 

e 

The solvability condition for the 0(e)-equation 

iwOX + lOxUir^+P {6YU(3n)Y + 

POyuoYri+P {dyf uinn + vui = , 

regarded as a differential equation for ui with respect to 
the variable 7] is 

Aqx + POyAoy + /3(0yAo)y = . (3) 

Adding Eq. (O, multiplied by (i/e)Ao exp((i/e)6'), to 
Eq. ([3]), multiphed by exp((i/e)0), we obtain 
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Uqx + ^POyUqy - -P{dY)'^Uo + POyyUq - -I^Uo = , 

e e 



2 



or, in initial coordinates {x, y) and introducing the wave 
number k = Oy , which is 0(1) in the method used, we 
get finally 

uqx + ^pkuoy — ifik^uo + (3kyUo — it^uo = . (4) 

The system of Eqs. ([2|) and jl]) describes the geometric 
optic approximation to the Eq. fl]), where two types of 
waves exist: propagating in positive direction along the 
y-axis, when A: > 0, and in negative direction when k < 
(the turning points are excluded from consideration in 
this paper). So we can use Eq. (|3]), replacing wo by u, 
as an approximate non-reflecting boundary condition at 
the boundaries of the form {{x,y)\y = const}. 

Formally we set the boundary value problem for Eq. ([T]) 
in the strip {x, y\x>0,a<y<b} (the initial- 
boundary value problem in the domain a < y < 6, if 
X is considered as the evolution variable), specifying the 
following boundary conditions: 

u — uj at a; = , 



Ux — 2(3\k\uy — ipk^u — I3\k\y — \vu = at y — 
Ux + 2l3\k\uy — ipk'^u + (3\k\y — ivu — Q at y — 



(5) 



where k = Oy for 6*, which is a solution of the Cauchy 
problem for the Hamilton- Jacobi equation ^ with the 
initial condition 

e{X,Y) = ei{Y) at x = o, 

specified for all values oiY . This initial condition appears 
as a result of the representation of initial data uj in the 
rapidly oscillating form 

ui = Ai{ey) exp((i/e)6l/(ey)) . 

Such a representation is not unique and depends on aux- 
iliary information that is used for recognizing the small 
parameter e, splitting the initial condition into the am- 
plitude and rapidly oscillating parts and determining the 
way of extrapolation of 9i to all values of Y and v outside 
the strip (this is needed to set the Cauchy problem for 
Eq. ^) 

Now we compare Eq. (j4|) with equations obtained by 
the formal factorization of Eq. ([1]) : 



vv 




ox oy I 



(6) 



which has an approximate character, because the opera- 
tors y/Tj3d/dy and ^Jdjdx^--\ii do not commute in gen- 
eral. The Pade approximant [3] Pq of the square root 
yjdjdx — \v with respect to djdx is 



dx 



1 + — — 

2v dx 



and using it in factors of Eq. ([6]), we obtain the equations 
Ux T 2,\f^Uy — \vu = . (7) 

If we put that the phase function does not depends on 
X, dx = 0, then from the Hamilton- Jacobi equation we 
set 



± 



and after substitution of these expressions into Eq. ^ we 
obtain Eq. ([7]) . Note that if the potential v is vanished at 
the boundaries, the boundary conditions ([7]) degenerate 
to the Dirichlet conditions. 

Because of that in Shibata's paper yj] was used an- 
other linear approximation of the square root, namely 
the linear interpolation between two points chosen with- 
out sufficient physical justification. 

Nevertheless, we will see later that our Eq. (|3|), having 
as an approximation the same linear nature, works quite 
well. 

Returning to the multiple-scale expansion, we consider 
the next step, which leads to the generalization of the 
rational-linear approximation discussed in Kuska's paper 

:5]. 

At O(e^) we obtain 

fidyu\y^+(i [Oyf U27)7) + VU2 + /SuoYY = , 

the solvability condition for which is 

iuix + /3 {0yui,^)y + I3dyuiyrj + fiu^yy = . (8) 

Putting, on the strength of the same argumentation as 
in the ug-case, 



ui = Aie'" , 

we get from Eq. ([5]) 

Aix+HGYAi)Y + POyAiy 



iPAoyy = . (9) 



In the same manner as Eq. ([4]) was derived, we obtain 
from this 



uix + 2[3kuiy—ij3k'^ ui + PkyUi — ivui — 

d[3Aoyy exp(- 6'(X, r)) = . 
e 



(10) 



We introduce now the approximation of the first order as 

U — Uq -\- Ui 

and obtain the equation for this quantity from 
Li/5(Eq. O) + eLi?S'(Eq. 

ehpAiyy exp(^ e{X, Y)) + 0{€^) = , 



3 



where LHS means 'left-hand side', as usual. Ex- 
pressing the terms i/3Aoyy exp(i 6'(X, F)) and 
\(3AiYY &^v{\S{X,Y)) from differentiated with re- 
spect to Y Eqs ([3|) and ([9]), after some algebra we 
obtain 



i(i/-h3/3fc2)My 



IVyU — (3kyyU " 3i(3kkyU — 6f3kyUy + 0{€'^) — . 

Then, the equations 
i{v + 3f3k'^)uy — Uxy T 3i|/c|Ma; =F |fc|(/3fc^ -I- 31^)^-1- 

iVyU ± l3\k\yyU it 31/3 | fc | | fc | u ±Q(3\k\yUy = 



(11) 



can be proposed as the corresponding non-reflecting 
boundary conditions, where the signs '— ' and cor- 
responds to y = a and y — h respectively. 

The rational-linear boundary conditions from Kuska's 
paper, written in our notations, read 

\{v + il3k^)uy-UxyT'iAk\u^T\k\{l3k^ + iv)u = 0. (12) 

These conditions were derived there by the factorization 
method with the Pade approximation !3!| Pi of the square 
root yjd/dx — \v. The last three terms of Eq. (fTT|) is 
absent in Eq. (|12p because fc is a constant in that paper. 
The term ivyU is absent due to the approximate character 
of the factorization method, mentioned above. 



III. NUMERICAL EXPERIMENTS 

As examples of application of the boundary conditions 
([5]) and (jlip we will present the numerical simulation of 
the Gaussian beam 



u(x,y) = 



XQ . y^ -yxo{2y+px) , 

■exp 1 j-^ — ; ^ , 13) 

X + Xq ' Al„ , „ \ I 



4:{x + xo) 



which was used also in the works 0, Si- As it is easily 
seen, the function from Eq. p3p is an exact solution of 
Eq. H]) with /3 = 1 and = 0. 

Choose xo = — ia, where a is a real number greater 
than zero. For the initial condition we have the expres- 
sion 



uq = u{0,y) = exp 



r 

'4a 



)=p(-.f) 



and we choose as the initial phase 



1 



% = 9{0,y) = --py. 



(14) 



The solution of the Cauchy problem for the Hamilton- 
Jacobi equation ^ with the initial condition is (see, 
e.g. i) 



1 



1 




FIG. 1: Contour of the logj^g(|ii|) for the boundary conditions 
(0). a = 2, p = 5, E/Eo = 1.216 . 10~*. 




FIG. 2: Contour of the logj^Q(|it|) for the boundary conditions 
(fTTI). a = 2, p = 5, E/Eo = 7.875 . 10~^ 



and this solution we use in the boundary conditions ([5]) 
and pTjl . Note also that in this case our first-order 
boundary conditions reduce to the Kuska's boundary 
conditions Eq. (jl2|l . The formal expansion of Eq. (fT3|) 
in powers of shows that in this case the latter can 

be considered as a small parameter, which is confirmed 
by the results of calculations. 

The calculations were done with the use of the Crank- 
Nicolson finite-difference scheme [l^ with the parameters 
of the Gaussian beam a ~ 2, p ~ 5 and a — 1/16, p — 40. 
The results of calculations for the first case with the grid 
size 1025 x 1025 are presented in FIG. [1] and FIG. H In 
FIG. [3] are presented the results of calculations with the 
use of the boundary conditions of Baskakov and Popov 
@, which are analytically exact, so this figure shows the 
effect of the used discretization. In the captions of these 
figures we present also the values of the relative energy 
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FIG. 3: Contour of the logjQ(jw|) for the Baskakov- Popov 
boundary conditions 0. a = 2, p = 5, E / Eq = 3.251 ■ 10"^. 

of the reflected waves 

PVmax /"Vmax 

E/Eq = / I y)\'dy/ / \u{0,y)\^dy. 

The results of calculations for the narrow Gaussian 
beam with parameters a = 1/16, p = 40, which were 
done in the domain [0 < x < 0.15] x [—2 < y < 2], show 
the following values of the relative energy: 

• 3.585 • 10-^ 5.677 • lO'^ and 1.834 • 10"^ for the ze- 
roth order, first order and Baskakov-Popov bound- 
ary conditions respectively on the grid 1025 x 1025; 

• 1.066 • IQ-'', 2.427 • lO^^ and 1.971 • lO^^ for the ze- 



roth order, first order and Baskakov-Popov bound- 
ary conditions respectively on the grid 513 x 513. 

These results confirm that the inverse width of the beam 
1/y^ plays the role of the small parameter and show 
also that for narrow beams the zeroth order boundary 
conditions can be better than the first order ones, even 
on the big grid, and more robust with respect to the 
roughness of the grid. 

The dependence of the accuracy of the boundary con- 
ditions (fT2l) on the beam width was investigated in some 
details in [5|. 

IV. CONCLUSION 

In this paper the absorbing boundary conditions ([5]) 
and (fTTI) for the variable coefficient Schrodinger type 
equation ([1]) were derived by the multiple-scale method. 
These boundary conditions explicitly take into account 
the variability of the coefficients and are easy to use. 
The solution of the Hamilton-Jacobi equation, which is 
required in these conditions, can be in most cases ob- 
tained analytically by the far-field approximations or nu- 
merically by the method of the eulerian geometric optics 

m- 

The reported method can be easily generalized to the 
many-dimensional case. 



Acknowledgments 

This work is supported by the Program No. 14 (part 
2) of the Presidium of the Russian Academy of Science. 



[1] Babich, V. M., and Buldyrev, V. S. (1972): Asymptotic 
Methods in Short- Wave Diffraction Problems. Nauka, 
Moscow [English translation: Springer Series on Wave 
Phenomena 4. Springer, Berlin et. al., 1991], Zbl. 
255.35002 

[2] Engquist, B. and Majda A. Mathematics of Computa- 
tion, 31, 639, (1977). 

[3] Bender, C. M. and Orszag S. A. Advanced Mathemati- 
cal Methods for Scientists and Engineers. McGraw-Hill, 
(1978). 

[4] Shibata, T. Physical Review B, 43, 6760, (1991). 

[5] Kuska, J.-P. Physical Review B, 46, 5000, (1992). 

[6] Baskakov, V.A. and Popov, A.V. Wave Motion, 14, 123, 



(1991). 

[7] Nayfeh, A.H. Perturbation methods. N.-Y.: John Wiley 

& Sons, (1973). 
[8] Trofimov, M.Yu. Technical Physics Letters, 31, 400, 

(2005). 

[9] Maslov, V. P. Russian mathamatical surveys, 42, 43, 
(1987). 

[10] Potter, D. Computational physics. N.-Y.: John Wiley & 
Sons, (1973). 

[11] Benamou, J.-D. Comm. Pure Appl. Math., 52, 1443, 
(1999). Also see INRIA report No. 4628, (2002). 



